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Abstract 

Wc study the isotropic Heisenberg chain with nearest and next-nearest neighbour 
interactions. The ground state phase diagram is constructed in dependence on the 
additonal interactions and an external magnetic field. The thermodynamics is stud- 
ied by use of finite sets of non-linear integral equations resulting from integrabiliy. 
The equations are solved numerically and analytically in suitable limiting cases. We 
find second and first order transition lines. The exponents of the low temperature 
asymptotics at the phase transitions are determined. 

1 Introduction 

Low dimensional quantum systems are of considerable current interest. On one hand they 
can be studied in experiments, where they are realised as quasi ID or 2D subsystems. On 
the other hand, some of the ID systems, like the spin-^ Heisenberg chain, can be solved 
exactly or otherwise non-perturbativly. The method of exactly solving quantum spin 
systems via Bethe ansatz is essentially restricted to ID models, but allows for solving 
models of coupled chains, see [13] and references therein. Depending on the topology, 
coupled chains (or spin ladders) may be considered as interpolations between ID and 
2D systems, or in the case of spin ladders with zigzag interactions the system may be 
viewed as a single chain with longer range interactions. 

In many cases these models are only studied in the ground state and without an 
external magnetic field, see e.g. [8]. Considering nonzero temperature and an external 
magnetic field is of interest for two reasons. First, the magnetic field can lead to (several) 
quantum phase transitions in the ground state, see e.g. [3]. Second, nonzero temperature 
and magnetic field are required for comparison with experimental work. 

The main goal of this paper is to study the thermodynamical properties of two 
quantum spin chains with competing interactions in an external magnetic field. Both 
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are generalisations of the standard spin- 2 Heisenberg chain and can be solved exactly 
via Bethe ansatz [14]. Here we investigate in more detail the thermodynamics which 
also leads to new information on the ground state. 

The paper is organized in the following way. First we introduce in sect. 2 the Hamil- 
tonians of the two models investigated in this paper. Then in sect. 3 we show the relation 
of these Hamiltonians to a row-to-row transfer matrix and present the definition of a 
quantum transfer matrix which allow for exactly solving the models via Bethe ansatz. In 
sect. 4 we derive non-linear integral equations determining the thermodynamical prop- 
erties of the models and consider their zero temperature limit. Using these equations 
we discuss the ground state phase diagrams in sect. 5 and present several results for the 
magnetic susceptibility, magnetisation and specific heat in sect. 6. Finally, we summarise 
our results. 



2 Hamiltonians 

In this paper we investigate the properties of two systems. First, the Hamiltonian of 
the Bethe ansatz solvable isotropic spin-^ chain with nearest neighbour interactions 
and three spin interactions between successive spins, referred to as next-nearest neigh- 
bour model or "model with NN interaction" [14, 2, 12], can be written in the form 
'Hnn = J'Hi + 02/^7^2 + T~ih, where 

i=i ^ ^ 

is the Hamiltonian of the standard spin-^ Heisenberg chain, and 

L 

n2 = Y,Si- X (2) 

i=l 

contains the three spin interactions. The Zeeman term 

L 

nh = -hY,st (3) 

takes account of an homogeneous external magnetic field h. The coupling 02 determines 
the relative strength of the three-spin interactions. 

The second system considered in this paper contains nearest and next-nearest neigh- 
bour interactions as well as four-spin interactions (in the following referred to as "model 
with NNN interactions" [14, 9]) which has the form 1~Lnnn = J'Hi + a^J^'Hz + 'Hh where 

^3=E 



i=l 



-SiSi+i + -SiSi+2 + 2 (Si^iSi+i) (SiSi+2) — 2 (SiSi+i) + - 



2 

(4) 

The operators Hi, {i = 1,2,3) and Tih commute mutually as well as with a transfer 
matrix t(A) constructed in the next section. These properties allow for exactly solving 
the models via Bethe ansatz. 
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3 Transfer matrices 



The i?-matrix belonging to the Heisenberg chain is given by 

^^J(A,^) = 5?5^ + (A-^)W (5) 



which is a solution of the Yang-Baxter equation. Here the indices in the first column 
denote states in the auxiliary space and the indices in the second column denote states 
in the quantum space. For a more detailed description of the notation see [1]. With 
)^(A,M) = i?;j(A,/x)e4 



(Lj)'g(A,^) = i?^J(A, /x)ej^ the operator defined on a chain of L sites 



t(A) =tr[LL(A,0)---Li(A,0)] (6) 

yields a family of commuting transfer matrices [t{X),t{fj,)] =0 for arbitrary X, fi € C The 
operators {TCi, i = 1,2,3) are given as logarithmic derivatives of this transfer matrix at 
the shift point A = 

Hi = ^r«(0)-§ (7) 
^2 = ^r(^^(0) + ^ (8) 
^3 = -^^^'^(0) + ^ (9) 

where r(A) = Int(A). 

Next, we introduce R matrices i?^J(A, fi) = -Rjf (a*, A) by a clockwise and i2^J(A, fi) = 
A) by an anticlockwise rotation and define a transfer matrix i(A) in the same way 
as t above. The partion function of the models without magnetic field can be expressed 
as 

N/2 



Z = lim tr 



i=l 



(10) 



with appropriate spectral parameters Uj [6], depending on the model. ^ The column- 
to-column transfer matrix of the corresponding two dimensional L x N lattice is called 
quantum transfer matrix (QTM). It is defined by 

t«™(A) = tr (z)L^™(A,0)L«™(A,n^/2) • • • L«™(A, 0)L«^^^(A, (H) 

where the magnetic field h is included by means of twisted boundary conditions via the 
diagonal matrix D = diag (exp(/3/i/2), exp(— /3/i/2)) and 

' '^^ l<^(A,^)e4 = i?5°,(M,A)e4, , odd. ^ ^ 



^Note that we ignored in (10) the additive constants of equations (7)-(9). 
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The monodromy matrix corresponding to the QTM is a representation of the Yang- 
Baxter algebra with intertwiner R. The partion function is given by 



ZN=tr{t^™{0)f . (13) 
Hence the free energy in the thermodynamic hmit is determined by 

f = -T hm lnA^™(0) (14) 

N—^oo 

where is the largest eigenvalue of the QTM. 

4 Non-linear integral equations 

For the standard spin-| Heisenberg chain one can derive different sets of non-linear 
integral equations (NLIE) determining the thermodynamical properties. Historically 
first, an infinite set of equations via TEA [11] was obtained. Then, as a second possibility, 
a set of only two equations [4, 5] was derived. In fact, it is also possible to find an arbitrary 
number of equations interpolating between these extreme schemes [10]. 

For the models investigated in this paper the set of two coupled non-linear integral 
equations was derived in [14] where also certain parameter ranges were treated numer- 
ically. Here, we performed further numerical studies of these equations and found that 
they are not valid for low temperatures in the vicinity of the phase coexistence if straight 
integration contours are used (see sect. 5). This is due to the fact that the imaginary 
parts of some Bethe ansatz numbers grow strongly, leading to a crossing of the integration 
contours by singularities of the integrands. 

To determine the free energy also for the cases where the two NLIE with standard 
contours are not valid, it is useful to utilise the fusion hierarchy of this model. One 
obtains in the usual way an infinite set of NLIE. The first equation is 

In yi {x) = -Z{x) + {s* In Y2) {x) (15) 

cosh(7rxj 

with * denoting convolution^ and six) = — - being the integration kernel. Z(x) 

2cosh(7rx) 

depends on the model and is given by 



Z{x) 



a2V^[5 2 ) fo'^ model with NN interactions , 

cosh {ttx) ,^r.\ 

.,2tanh2(7rx)-l , , . ^ ' 

a^v p — — r , for the model with NNN mteractions . 

cosh(7rx) 

The other equations are independent of the model and read 

lny,(x) = (s * In {Yj.^Yj+i)) {x) . (17) 



\9*h) {t)=Sr^g{t-T)h{r)dr 
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The magnetic field does not enter explicitly in these equations, it only fixes the asymp- 
totic behaviour of the y-functions. For zero magnetic field it reads ii^\x\^ooyji^) = 
j{j + 2) and for /i / 

/^i+i _ 2-0+1) \^ 
lim yj{x) = - 1 (18) 

x|— »oo \ Z — Z I 

with z = e^'*/^. Note that the asymptotic behaviour and therefore the system of equa- 
tions is invariant under a change of sign of the magnetic field. The free energy is given 
by 

oo 

/ = e„-r/dx-i!i4M^ (19) 

J 2cosh(7rxj 

— oo 

with Co = — J In 2 for the model with NN interactions and cq = — J In 2 + a3^^C(3) for 
the model with NNN interactions, C denoting the Riemann (^-function. 

It is possible to close the set of infinitely many integral equations after the {k — l)th 
equation. The minimal paramter k = 1 gives the set of two NLIE presented in [14]. Next 
we show the results for k >2 following [10]. 

One can find suitable functions b, b, !B(x) := l-|-b(x) and *B(x) := l + b(x) satisfying 
!B(x)!B(x) = Yii.{x) and therefore yielding the functional relation 

yfc„i(x + ^)yfc„i(x - ^) = yfc_2(x)«(x)»(x) . (20) 

For (20) the cases k = l and k = 2 are exceptional. For k = l the equation is not used, 
for k = 2 we use Yq = 1. This leads to the following NLIE 

lnyfc_i(x) =s*ln (^^2) (x) + s * In (QS^) (x) . (21) 

Finally, two NLIE equations for b, b close the equation system exactly 

Inb(x) = ^ + (s *lnyfe_i) (x) + (K*ln<B) (x) - (K*ln^) {x + i) (22) 

_ Bh 

Inb(x) = --^ + {s*\uYk-i){x) + {K*\u^){x)-{K*\u^){x-i) . (23) 

In (21) for k = 2 and in (22), (23) for A; = l, instead of s *lnlo the inhomogeneity of (15) 
has to be used. The kernel function k is given by 

k{x) = — / dk————e'^''. (24) 
^ ' 2tt J 2cosh(A:/2) ^ ' 

—00 

Here, the magnetic field enters the equations explicitly. The different sets of equations 
for different k are all equivalent. The equations for larger k allow for the investigation 
of lower temperatures even with straight contours. 
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Calculating the zero temperature limit of the fusion hierarchy with = Tlnyj in 
the usual way [11] one obtains only one dressed energy e := £i which can have negative 
values, all other dressed energies are strictly positive and hence do not contribute. The 
dressed energy e is determined by 

e(A) + ljdfi ^ ^ eip) = eo(A) (25) 



with 



and = {/i € ]R|e(/i) < 0}. Note that only the absolute value of the magnetic field 
enters the equation (25) via the bare energy (26). This is due to the fact that only the 
absolute value of h enters the fusion hierarchy via the asymptotic behaviour (18). This 
statement becomes even more obvious if one derives (25) from the system of only two 
NLIE. Here one of the steps of the derivation uses the observation that the function 
b drops out for > and b drops out for /i < and in both cases (25) is obtained. 
Furthermore a rescaling of the couplings 02, as was applied by a = a2V for the model 
with NN interactions and a = a^v"^ for the model with NNN interactions. This will be 
useful in the discussion of the ground state phase diagram in the next section as in the 
new variable certain critical points occur at a = ±1. 

The dressed energy can be solved in two limiting cases. The first one is at h = and 
\a\ < 1 by Fourier transform 

£(X) = l~ ^^4^ + " tanh(7rA)] for Hnn , ^^7) 

1 ~ cosh(7rA) [1 - " + 2a tanh^(7rA)] for Hnnn ■ 

The second case is for h > hf where the integral vanishes and e(A) = £o{X). hf is the 
saturation field, i.e. the value of the magnetic field corresponding to the phase transition 
into the ferromagnetic phase. 



5 Phase diagrams of the ground state 

In this section we first give our results on the ground state phase diagram for both models 
(Fig. 1). Our results differ a little from those of [14]. For all couplings a there is a phase 
transition into a ferromagnetic phase and phase coexistence for the lines \a\ > 1, h = 0. 
But only in the model with NNN interactions and positive coupling a we find a phase 
transition between a commensurate and an incommensurate phase. 

Looking at the dressed energy for the models with zero magnetic field one can under- 
stand the phase diagrams. The phase transitions in dependence on the magnetic field h 
correspond to the opening and closing of Fermi seas, i.e. appearance or disappearance of 
intervals of negative energy modes of the dressed energy. Qualitatively, these transitions 
occur at magnetic fields h coinciding with extremal values of the dressed energy at zero 
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(a) Model with NN interactions (b) Model with NNN interactions 



Figure 1: Phase diagrams of the model with NN and NNN interactions. Both phase 
diagrams are symmetric with respect to a change of sign of the magnetic field h, and for 
the model with NN interactions (a) also with respect to a change of sign of the coupling 
a. In these pictures the normalization J = 2 in Hnn and TCnnn is used. The first order 
phase transitions are exactly at h = 0. 

field as plotted in Fig. 2. This is very much like the discussion of van-Hove singularities 
of free particle systems. However, here we deal with an interacting system for which the 
chemical potential is not identical to h but equal to \h\. So only the negative extremal 
values of the dressed energy are relevant. Although the analytic solution (27) is strictly 
valid only for \a\ < 1 we use these formulas for slightly larger values of |a| where they 
should be good approximations to the true solution. 

Depending on the longer range coupling a, the model with NNN interactions has 
one or two Dirac seas. For a < 0.2 only a second order phase transition into the 
ferromagnetic phase and a first order transition line with non vanishing spontaneous 
magnetisation for a < —1, h = exist (Fig. 2c and 2d). For a > 0.2 (Fig. 2b) the dressed 
energy possesses two local minima and one local maximum. The two local minima have 
identical value as e(A) is an even function. So here a phase transition between two 
phases exists (denoted "commensurate" and "incommensurate") . The tricitical point can 
be determined to a = ^ ~ 0.2056 and hf = ^ hy analysing the solution of the dressed 
energy in the ferromagnetic phase. 

The phase diagram of the model with NN interactions is symmetric under a sign 
change of the coupling a, as in the NLIE this can be compensated for by a sign change 
of the spectral parameter A. Hence for the model with NN interactions it is sufficient 
to look at couplings a > 0. This model has always only one Dirac sea (Fig. 2a), hence 
there is only a second order phase transition between the antiferromagnetic and the 
ferromagnetic phase and a first order transition line with non vanishing spontaneous 
magnetisation for \a\ > 1, h = 0. 

The critical field corresponding to the phase transition from the antiferromagnetic 
into the ferromagnetic phase can be determined as usual [7] from the bare energy (26). 
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(a) NN model with positive a (b) NNN model with positive a 




(c) NNN model with negative a (d) enlarged case a = —1.1 of (c) 



Figure 2: The dressed energies according to equation (27) for h = and normalisation 
J = 2. The data are approximate for \a\ > 1. 



For the model with NNN interactions one obtains 
'2J-ai|j, 

2J 



h 



f 



+ 16J- 



3,/24- 



24^+4-4-72^ 



24- 



24-%-+4 -24- 



24- 



24^+4 -24A 



for a < 7rV48 , 
for a > TT^/AS. 



(28) 



For a < this is a linear relation between the critical field hf and the coupling a. 

For a > 7r^/48 the relation is non-linear, however with linear asymptotic behaviour for 
large values of a 

(29) 



h 



f 



J 



a 



The applicability of (28) is restricted to really large values of q, e.g. the error gets smaller 
than 1% for a ~ 5.8. For the model with NN interactions the equation determining the 
minimum of the bare energy is cubic, whereas it is biquadratic for the model with NNN 
interactions. For this reason we want to give only numerical values for the critical field 
hf for the model with NN interactions. These can be taken from Fig. la. 
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(a) x{h) for a 



-1, r = 0.01 



(b) x(T) for a 



-1, /i ~ 



Figure 3: The magnetic susceptibility for the model with NNN interactions and negative 
coupling a ~ — 1. 



6 Magnetic susceptibility, magnetisation and specific heat 

In this section we present the magnetic susceptibility x foi' typical values of the coupling 
a in dependence on the magnetic field and temperature. We also show evidence that 
the phase transition at /i = is of first order. Finally, we calculate the specific heat c at 
the phase transitions. The derivatives of the free energy are obtained by differentiating 
(19) and deriving new integral equations for the logarithmic derivatives of the auxiliary 
functions involving the auxiliary functions as external parameters. 

As the situation in the model with NN interactions is very similar to the one for the 
model with NNN interactions and negative coupling a, we will focus in the following on 
the model with NNN interactions and only sometimes give comments on the model with 
NN interactions. For the numerical calculations we always use the normalisation J = 2. 

For the model with NNN interactions the magnetic susceptibility is shown for typical 
values of a. In Fig. 3a x(^) is shown for a = —l and r = 0.01. One sees that there is a 
maximum at /i ~ 7.2 corresponding to the phase transition between the antiferromag- 
netic and ferromagnetic phase and another maximum at /i = also corresponding to a 
divergence at r = (Fig. 3b). This qualitative picture is also true for a < —1. 

For larger values of a the low field maximum exists until a<— 0.2 and h <^ hf but it 
does not correspond to a phase transition because the magnetic susceptibility does not 
diverge for T — > 0. 

For a > 0.206, the value corresponding to the tricritcal point, two maxima exist 
(Fig. 4a). Here also the maximum at lower magnetic field diverges for T ^ and 
so corresponds to a phase transition (Fig. 4b). The value of the lower critical field 
decreases with increasing a and turns at a = l. This fact, known from the functional 
behaviour of the dressed energy is also supported by the numerical data for the magnetic 
susceptibility at finite temperature. In Fig. 4c one clearly sees that the maximum of the 
magnetic susceptibility occurs at a finite magnetic field for a < 1. 

The value of the magnetic field corresponding to the phase transition between the 
commensurate and incommensurate phase could not be determined analytically. How- 
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(b) x{T) for a = 0.7 



a = -1 - 

a = -1.01 

lc-04 0.001 
h 

(c) xW for T = le-3 (d) m{h) for T=le-5 

Figure 4: The magnetic susceptibility for the model with NNN interactions in figures 
(a)-(c) and the magnetisation in figure (d). Fig. (b) clearly shows that for a = 0.7 the 
magnetic susceptibility diverges for /i = 1.455, which is in contrast to the cases with lower 
and higher magnetic fields. 

ever, good numerical results are obtained by calculating the magnetic susceptibility at 
finite but low temperature for different values of the magnetic field and determining the 
local maximum. Doing this for different temperatures one can also estimate the error. 
The line between the commensurate and incommensurate phase in Fig. lb is located in 
this way, where errorbars are within linewidth. 

The magnetisation for /i — >■ is shown for the model with NNN interactions and 
a ~ — 1 in Fig. 4d. Clearly for a < — 1 and low temperature the magnetisation has a 
finite asymptotic limit for small magnetic fields, whereas this is not the case for a = — l. 
So, a first order phase transition exists for a < — 1. This statement also holds for h = 
and a > 1 and for the model with NN interactions with zero magnetic field and \a\ > 1. 

Finally, we determine the asymptotic behaviour of the specific heat for T — > at the 
phase transitions. We find the specific heat for h = and \a\ = 1 to very low temperatures 
and with very high accuracy. 

The specific heat vanishes as T^/^ for the model with NNN interactions at a = — 1 
(Fig. 5a) and as T^/^ (Fig. 5b) for a = 1. For the phase transition into the ferromagnetic 
phase and with slightly lower numerical accuracy for the phase transition between the 



(a) x{h) for a = 0.7, T = 0.01 
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(c) Q = ^; ft = f ; X = 0.2499 ± 5 (d) NN-model: a = 1- h = Q; x = 0.33339 ± 1 



Figure 5: The specific heat for the model with NNN interactions (Fig. (a)-(c)) and the 
model with NN interactions (Fig. (d)). The low temperature exponent x of T (c ~ T^) 
is determined by a fit on the numerical data. The error refers to the last digit and only 
includes the statistical error of the fit. 

commensurate and incommensurate phase we find the asymptotic behaviour T^/^. At 
the tricritical point the specific heat vanishes like T^l^ (see Fig. 5c) as predicted in [14] . 
For the first order phase transitions we were not able to compute the specific heat at 
sufficiently low temperatures with sufficiently high accuracy to find consistent results for 
the low temperature asymptotics for |a| > 1. Here further investigations are necessary. 

Again the situation in the model with NN interactions corresponds to the one in the 
model with NNN interactions and negative coupling a. In particular this means that 
the specific heat at a = 1 and /i = vanishes like T^/^ which is shown in Fig. 5d. 

7 Conclusion 

We studied the thermodynamics and ground state phase diagrams of two integrable 
models containing the standard spin-^ Heisenberg Hamiltonian and additional competing 
interactions. 

The ground state phase diagrams depending on the external magnetic field and the 
longer range coupling a are constructed. They contain ferromagnetic and antiferromag- 
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netic phases. In both models there exist second order phase transitions between these 
phases and first order phase transition hnes with non- vanishing spontaneous magnetisa- 
tion. Only the model with NNN interactions with positive a contains a phase transition 
between a commensurate and an incommensurate phase. 

The NLIE describing the models at finite temperature are solved numerically for 
typical values of the coupling a and the magnetic field h. 

The vicinity of the phase coexistence in both models is difficult to investigate. For 
sufficiently low temperatures the NLIE are numerically ill-posed if straight integration 
contours are used. For reaching low temperatures, either the contours have to be de- 
formed - or as chosen in our approach - the truncation level has to be increased. 



Acknowlegment 

The authors like to acknowledge support by the research program of the Graduiertenkol- 
leg 1052 funded by the Deutsche Forschungsgemeinschaft. 



References 

[1] F.H.L. Essler, H. Frahm, F. Gohmann, A. Kliimper, and V.E. Korepin, The one- 
dimensional Hubbard model, Cambridge, 2005. 

[2] H. Frahm, Integrahle spin-1/2 XXZ Heisenberg chain with competing interactions, 
J. Phys. A: Math. Gen. 25 (1992), 1417-1427. 

[3] H. Frahm and C. Rodenbeck, A generalized spin ladder in a magnetic field, 
Eur. Phys. J. B 10 (1999), 409-414. 

[4] A. Kliimper, Free energy and correlation lengths of quantum chains related to re- 
stricted solid-on-solid lattice models, Ann. Phys. 1 (1992), 540-553. 

[5] , Thermodynamics of the anisotropic spin-1/2 Heisenberg chain and related 

quantum chains, Z. Phys. B 91 (1993), 507-519. 

[6] A. Kliimper and K. Sakai, The thermal conductivity of the spin-1/2 XXZ chain at 
arbitrary temperature, J. Phys. A: Math. Gen. 35 (2002), 2173-2182. 

[7] V.E. Korepin, N.M. Bogoliubov, and A.G. Izergin, Quantum Inverse Scattering 
Method and Correlation Functions, Cambridge, 1997. 

[8] M. Miiller, T. Vekua, and H.J. Mikeska, Perturbation theories for the s = spin 
ladder with a four- spin ring exchange, Phys. Rev. B 66 (2002), no. 13, 134423. 

[9] N. Muramoto and M. Takahashi, Integrable Magnetic Model of Two Chains Coupled 
by Four-Body Interactions, J. Phys. Soc. Jpn. 68 (1999), no. 6, 2098-2104. 

[10] J. Suzuki, Spinons in magnetic chains of arbitrary spins at finite temperatures, 
J. Phys. A: Math. Gen. 32 (1999), 2341-2359. 



12 



[11] M. Takahashi, One-Dimensional Heisenherg Model at Finite Temperature^ 
Prog. Theor. Phys. 46 (1971), no. 2, 401-415. 

[12] A.M. Tsvelik, Incommensurate phases of quantum one- dimensional magnetics, 
Phys. Rev. B 42 (1990), no. 1, 779-785. 

[13] A. A. Zvyagin, Bethe ansatz solvable multi-chain quantum systems, J. Phys. A: 
Math. Gen. 34 (2001), R21-R53. 

[14] A. A. Zvyagin and A. Kliimper, Quantum phase transitions and thermodynamics of 
quantum antiferromagnets with next-nearest-neighbor couplings, Phys. Rev. B 68 
(2003), 144426. 



13 



